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■ Abstract. We review some of the recent advances in nonlinear pulsa- 

Q\ ', tion theory, but also insist on some of the major extant shortcomings. 



1. Introduction 

lO 

■ It is a great pleasure for me to dedicate this lecture to Art Cox whose contribu- 

tions to nonlinear pulsations go back to the very pioneering years of numerical 
hydrodynamics some 35 years ago. It is fair to say that there isn't a type of star 
that Art hasn't attempted to model. 
| It is of course impossible to cover in this short review all the topics of 

fNj ■ interest to nonlinear stellar pulsations and I will have to make a selection that 

necessary reflects my biases. There has been a huge amount of work done on 
stellar pulsations and I also apologize upfront for many important omissions. 
Extensive references are provided in the excellent reviews of Gautschy and Saio 
(1995, hereafter GS). 

There are basically two approaches to nonlinear stellar pulsations that are 
complementary in some ways. The first, and oldest, is numerical hydrodynamics. 
While this is a 'brute force' approach, it has the advantage that with state-of- 
£^ ■ the-art physical input and numerical methods it can yield accurate information 

! about the nonlinear pulsations of individual stellar models. The second approach 

is the amplitude equation formalism, and it is perhaps of a more fundamental 
nature (e.g. Buchler 1988). It gives a broad overview of the possible behavior, 
^ ' such as modal selection (bifurcations in modern language) and the effects of 

resonances (Buchler 1993). We note that currently this is the only tool with 
which we can understand nonlinear nonradial pulsations. 

About 25 years ago John Cox (1975) felt that "overall, pulsation theory 
and its applications are in a fairly satisfactory state, except for a few disturbing 
problems" . In the intervening time, of course, a good deal of progress has been 
made. However, there remain some "disturbing problems", and in fact some 
new ones have appeared recently. 

Perhaps the most important progress came from outside pulsation theory, 
namely from a revision of the stellar opacities (Iglesias & Rogers 1992, Seaton, 
Kwan, Mihalas & Pradhan 1994), and it essentially solved two longstanding 
Cepheid problems. First, the so-called Cepheid bump mass problem (e.g. Art 
Cox 1980) essentially disappeared and the agreement with observation appears 
now to be quite satisfactory for the Galactic Cepheids. (Moskalik, Buchler & 
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Marom 1992, Kanbur & Simon 1993). Second, the beta cephei models finally 
became linearly unstable (cf. GS). 

The very first numerical Lagrangean hydrodynamical computations of Ceph- 
eid variables were already apologetic about the poor resolution of the partial 
hydrogen ionization front during the pulsation. In the late 70s Castor, Davis 
& Davison (1977) developed the first code that could track the moving sharp 
features and Aikawa & Simon (1983) started to use it systematically. More 
recently, taking advantage of new developments in computational physics and 
of faster computers, several groups have developed more flexible adaptive codes 
(Gehmeyr & Winkler 1992, Dorfi & Feuchtinger 1991, Buchler, Kollath & Marom 
1997). With these codes it is now possible to obtain a good spatial resolution 
that satisfactorily resolves all shocks and ionization fronts and achieves a much 
enhanced numerical accuracy of the pulsation. The most striking improvement 
is in the smoothness of the lightcurves and radial velocity curves. Fortunately, 
though, we do have to discard the Lagrangean results since quantities such as 
Fourier decomposition parameters are not substantially different from those ob- 
tained with Lagrangean codes. 

However, instead of congratulating ourselves on these and other successes, 
it is perhaps more useful to dwell on the " disturbing problems" . 

2. 'Disturbing Problems' 

2.1. Low metallicity Cepheids: 

The microlensing projects have provided us with a large treasure trove on vari- 
able stars in the Magellanic Clouds, and since these galaxies have been found 
to be metal deficient compared to the Galaxy the new observations have consid- 
erably enlarged our data base, especially for the Cepheids since they should be 
strongly affected by metallicity content. 

The Fourier decomposition parameters of the fundamental Cepheid light- 
curves in the SMC (Beaulieu & Sasselov 1997) and the LMC (Beaulieu et al. 
1995, Welch et al. 1995) indicate that the fai phase progression is very similar 
to that of the Galaxy, although it may be shifted by ± 1 or perhaps 2 days 
in period. However, the size of the excursion in <j>2\ in the resonance region is 
essentially the same as in the Galaxy. 

A comparison of the linear bump Cepheid models with mass-luminosity 
relations derived from evolutionary computations show up an irreconcilable dif- 
ference with the observations (Buchler, Kollath, Beaulieu & Goupil 1996). Fur- 
thermore the nonlinear calculation of low Z Cepheid model pulsations give <j>2\ 
in which the size of the excursion in the 10 day resonance region almost vanishes 
as one goes to Z values of 0.005 (Fig. 1) (see also poster by Goupil). 

As far as beat Cepheids are concerned, the hope had been that the observed 
period ratios could give a powerful constraint on the stellar model parameters. 
While globally, it might appear that the observed period ratios of the Galactic, 
LMC and SMC beat Cepheids are in agreement with the linear models obtained 
with the new opacities, when looked at in detail, i.e. by considering individual 
stars, the agreement is no longer as good. More seriously, (unknown) nonlinear 
period shift corrections of as little as 0.1% can give substantially different or 
uncertain mass assignments (Buchler et al. 1996). 
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Figure 1. Behavior of the Fourier phase 02i as a function of pulsation 
period for various values of metallicity. 

2.2. Overtone Cepheids (s Cepheids) 

The Fourier decomposition parameter cp2i for the first overtone Cepheids display 
a "Z" shape around a period of 3 days for the Galaxy (e.g. Antonello et al. 1990) 
and around ~ 3.5 and 4.0 days for the LMC (Beaulieu et al. 1995, Welch et 
al. 1997) and SMC (Beaulieu & Sasselov 1997), respectively. Hydrodynamical 
models of overtone Cepheids do not agree with the observations (Antonello & 
Aikawa 1993). Similar results were obtained by Schaller & Buchler in a fairly 
extensive survey of s Cepheids (unpublished preprint, 1994). 

2.3. RR Lyrae: 

Overall, the modelling of RR Lyrae pulsation gives decent agreement (except 
for 'double-mode' RRd pulsations) (cf GS for references), but when a more 
detailed comparison with observation is carried out in a systematic fashion, 
serious discrepancies pop up as Kovacs & Kanbur (1997) show. 

2.4. Population II Cepheids 

BL Her stars: The modelling of these low period stars (Hodson, Cox & King 
1982, Buchler & Buchler 1994; cf. also GS) shows overall agreement for the 
lightcurve data, but the 02i are considerably smaller than the observational 
data indicate. 

W Vir and RV Tau stars: It is now well known that the W Vir and RV Tau stars 
belong to the same group (Wallerstein & Cox 1984) and that the properties vary 
gradually from the low period, low luminosity W Vir stars to the long period, 
high luminosity RV Tau stars. Although the observations are not very extensive 
they indicate that the W Vir stars are periodic up to « 15 days from whence they 
start showing alternations in the cycles, alternations that become increasingly 
irregular with 'period'. The mechanism for this irregular behavior remained a 
mystery until relatively recently. 

Numerical hydrodynamical modelling of sequences of W Vir models (Buch- 
ler & Kovacs 1987) uncovered very characteristic nonlinear behavior that goes 
under the name of low dimensional chaos. Since the concept of low dimensional 
chaos is still very new in Astronomy we stress that this behavior is very different 
from a static multi-periodic, and also different from an evolving multiperiodic 
system. (For the reader familiar with chaos we mention that the presence of 
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period doubling along sequences of models quite clearly shows the presence of 
a horseshoe dynamics with almost one-dimensional return maps. The chaos in 
these models is of the stretch-and-fold type, very similar to the one that occurs 
in the Rossler system of 3 ODEs, e.g. Thompson & Stewart (1986). 

There is however a discordance with observations, that is the onset of period 
doublings and chaos occurs already at 7-10 days, rather than at the ps 15 days 
indicated by the observations. 

The numerical modelling of RV Tau behavior is much harder because the 
ratio of growth-rate/pulsation frequency is much greater. The pulsations are 
thus much more violent and result in sudden loss of the whole envelope in our 
calculations (see however Fadeyev & Fokin 1985, Takeuti & Tanaka 1995). We 
think that a physical dissipation mechanism is missing from our modelling, even 
though we solve the radiation hydrodynamics equations. Most likely turbulent 
dissipation plays a role in taming the pulsations of these stellar models. 

The occurrence of chaos in hydro dynamical models is quite robust as we 
have already indicated, but could it be an artifact of the theoretical modelling, 
even though it was confirmed with a totally different code (Aikawa 1990). Clearly 
it needed to be challenged by observation. A recent nonlinear analysis that goes 
under the name of 'global flow reconstruction' rather conclusively shows that 
the irregular lightcurve of R Set, a star of the RV Tau type, is the result of a 
low-dimensional chaotic dynamics (For details we refer the reader to Buchler et 
al. 1995, or to a didactic review Buchler 1997). More specifically, the analysis 
establishes that the dimension is as low as 4. In other words, the lightcurve is 
generated by 4 coupled ordinary differential equations. Put differently, if s(t) 
denotes the magnitude of the star, then at any time t n the lightcurve is a func- 
tion of only four preceding times, s(t n ) = F[s(t n -i), s(t n - 2 ), s(t n - 3 ), s(£ n _ 4 )] 
This result is quite remarkable since the pulsations of this star are quite vio- 
lent (factors of 40 changes in luminosity!) with shocks and ionization fronts 
running about. As a physicist, though, we are not satisfied merely with this 
result, but would like to know what more physics we can learn about this star. 
A four dimensional dynamics indicates that probably two vibrational modes are 
involved in the dynamics. This is strongly corroborated by a linearization of 
the dynamics about the equilibrium that tells us that two spiral stability roots 
are involved, one unstable with frequency /o=0.0068 d _1 , the other stable with 
frequency /i=0.014 d -1 ^ 2/o The physical picture that emerges is then that 
the irregular lightcurve of R Set is the result of the nonlinear interaction be- 
tween an unstable, lower frequency mode and a linearly stable overtone with 
approximately twice the frequency. 

A recent nonlinear analysis of the AAVSO lightcurve of AC Her similarly 
shows low dimensional chaos (Kollath et al. 1997). 

To summarize, the predictions of the nonlinear hydrodynamics, viz. that 
the irregular behavior of these stars is due to low dimensional chaos, are thus 
confirmed, but better numerical modelling is necessary to achieve closer agree- 
ment with observations. 

Could there be a common cause for most of the "disturbing problems"? It 
is very unlikely that the opacities or the equation of state are still at fault. A 
better treatment of radiative transfer (e.g. poster by C. G. Davis) is also not 
likely to fix most of the discrepancies. However, as we have already pointed 
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Figure 2. Behavior of the lightcurve as a function of zone number 
with the additional zones between 11,000 and 65,000 K; top: radiative 
adaptive code; bottom: turbulent diffusion code. 

out, we seem to be missing a physical dissipation mechanism in our radiative 
codes. In Lagrangean codes it is necessary to include pseudo-viscosity (a la 
von Neumann-Richtmyer) in order to handle shocks. While this viscosity is 
ok for many explosive shock problems, such as supernova explosions, it unfor- 
tunately also provides artificial and unphysical dissipation elsewhere. Kovacs 
(1990) found that when he reduced the artificial dissipation to a minimum the 
nonlinear pulsation amplitudes kept increasing to unrealistically large values. 
Similarly, when we increase the spatial resolution of the models the pulsation 
amplitudes also increase (Fig. 2). In fact it turns out that no combination of 
linear and quadratic viscosity parameters give satisfactory fundamental and first 
overtone models. (We hasten to add though that the Fourier decomposition pa- 
rameters and the Floquet stability coefficients of the limit cycles, fortunately, 
are reasonably independent of these changes, so that we do not have to throw 
away everything we have done so far!) 

3. Turbulent diffusion and convection 

I had always hoped that, at least near the blue edge of the instability, the major 
effect of convection was static and would thus merely cause a small systematic 
change in the structure of the models, so that we could get away with purely 
radiative hydro models. (Of course it is the important dynamic effect of con- 
vection which gives rise to the red edge.) The problems and tests that we have 
described above however indicate that we have to include turbulent convection 
in the hydrocodes in order to provide a missing powerful dissipation mechanism. 

Turbulent convection is of course a 3D phenomenon and at present, and 
for some time to come, it is not possible to run realistic 3D pulsation models. 
Progress has been made with relatively idealized 3D convection modelling, but 
it is slow and these calculations do not yet allow us to extract the ID recipes 
that we need for our radial pulsation codes. In the meantime we have to rely on 
ad hoc ID recipes with ad hoc parameters. 

The earliest models for convection were local both in time and in space 
and were found to be inadequate for stellar pulsations. Today we have a family 
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of time-dependent turbulent diffusion models that go back to Spiegel (1963), 
Unno (1967) and Castor (1968). A simplified version was implemented in a 
hydrodynamics code by Stellingwerf (1982). Recent applications have been made 
by several groups, viz. Bono et al. (1997), Gehmeyr (1992), Feuchtinger & Dorfi 
(1996) and by Yecko, Kollath k, Buchler (1997). The strategy has been and 
remains to compare the predictions of the models with observations and from 
thence calibrate the unknown parameters. 

Our own numerical testing shows that with a turbulent diffusion model 
the saturation amplitude of the pulsations becomes largely independent of the 
zoning (Fig. 2, kindly prepared by Phil Yecko). Another positive point is that 
with the reduced pulsation amplitudes the shocks are absent or much weaker. 
This in turn allows one to reduce the artificial dissipation to a very small value. 

As a word of caution, we note though that these models may still be too 
local in space (they only have a diffusion operator for the turbulent energy) 
and the existence of plumes may have to be taken into account as suggested by 
Rieutord & Zahn (1995). 

4. Amplitude Equations 

We have already pointed out that the amplitude equation formalism offers an 
alternative to 'brute force' numerical modelling. We would like to stress here 
that contrary to the claim of GS this formalism is not an Ansatz, but is a 
mathematically rigorous aproach, namely normal form theory. Essentially, the 
only restriction is that the formalism applies to weakly nonadiabatic pulsators. 
Many of the interesting stars, viz. the classical Cepheids, the RR Lyrae, the delta 
Scuti and the white dwarfs definitely fall into that category. For details of the 
formalism as applied to stellar pulsations we refer the reader to reviews (Buchler 
1988, 1993). We note also that the formalism has recently been extended to 
nonradial pulsations, in an Eulerian formulation by Goupil & Buchler (1994) 
and in a Lagrangean one by van Hoolst (1994). 

In a nutshell, the formalism reduces the PDEs of hydrodynamics and radi- 
ation transfer to a small set of ODEs for the amplitudes of the excited modes. 
The structure of the equations is uniquely determined by the types of resonances 
that occur among the linear modes of oscillation; the remaining physics is all 
contained in the values of the nonlinear coefficients. The amplitude equations 
are generic and capture the essence of the behavior of the system; it is not as- 
tonishing thus that they pop up in many different areas of physics, chemistry, 
biology, etc.. 

The solutions of the amplitude equations (usually the fixed points) tell us 
then about the possible types of behavior for a sequence or array of stellar 
models. They also explain the effect of resonances on the morphology of the 
lightcurves and radial velocity curves. Perhaps the most useful and best known 
application of the formalism has been to describe the behavior of the Fourier 
decomposition parameters through the Hertzsprung progression of the bump 
Cepheids (e.g. Buchler 1993). 

5. Potpourri 
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Figure 3. Schematic behavior of the fundamental and first over- 
tone amplitudes for nonresonant scenario; left: fundamental amplitude; 
right: first overtone amplitude. 
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Figure 4. Behavior of the fundamental and first overtone ampli- 
tudes for the fundamental and first overtone amplitudes of the SMC 
Cepheids; left: fundamental amplitude; right: first overtone amplitude 
(courtesy of Beaulieu). 



5.1. 'Double-mode' behavior 

As an application of the formalism to beat (double-mode) behavior we have 
plotted in Fig. 3 the predictions of the amplitude equation formalism in which it 
is assumed that two nonresonant modes (e.g. the fundamental and first overtone) 
interact and can give rise to a double-mode pulsation along a single parameter 
sequence of models (scenario 'AB' in Fig. 1 of Buchler & Kovacs 1986). Here 
we denote by Aq both the Fourier amplitude of the fundamental mode for the 
fundamental pulsators, and also the amplitude of the fundamental component in 
the case of double mode pulsators. The right figure shows the corresponding first 
overtone amplitudes. Note that the transition from single mode to double-mode 
occurs smoothly for both amplitudes. In a realistic sample of models one would 
of course get a dispersion both vertically and horizontally, but the conclusion 
remains unaffected. 

Let us compare this now to the observations, first to the RR Lyrae in M15 
of Sandage et al. as reproduced in Buchler & Kovacs 1986. In their fig. 8 the 
first overtone amplitudes are displayed on the left as dots for the RRc and as 
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Figure 5. Variation of the Fourier decomposition parameters over 
the Blazhko cycle for RR Lyrae and AR Herculis. 

crosses for the RRd, whereas the fundamental amplitudes are the dots on the 
right for RRab and open circles for the RRd. The A\ amplitudes of the RRc 
and RRd stars indeed vary continuously, but the Aq amplitudes of the RRd are 
considerably smaller than those of the RRab. We are forced to conclude that a 
nonresonant scenario is not in agreement with the observations. To explain the 
jump in the fundamental amplitudes from the RRab to the RRd is is necessary 
to invoke the presence of a resonance. (A jump in the amplitudes might also be 
brought about by higher order, viz. quintic nonlinearities, but in no studies so 
far have such nonlinearities ever been found to play a role, and furthermore the 
coefficient of the cubic nonlinearity would have to have a sign opposite to its 
usual one.) However, no low order resonances are present in the stellar parameter 
range of the RR Lyrae, and it therefore has to be a higher order resonance that 
is involved. We note in passing that therefore these stars were better called beat 
RR Lyrae because more than 2 modes are involved. 

Let us now turn to the Cepheids. J.-P. Beaulieu has kindly provided me 
with his SMC Fourier decomposition data that are displayed in Fig. 4. The 
first overtone amplitudes of the beat Cepheids fall right into the range of the 
s Cepheids, but the fundamental amplitudes again are much smaller for the 
beats than for the fundamental Cepheids. We are forced to interpret this to 
mean that a resonance must also be involved in the beat Cepheids. 

5.2. The Blazhko effect 

Several mechanisms for the Blazhko effect have been proposed (cf. GS), but a 
fully satisfactory understanding has so far defied us. In the following we want 
to present an observational constraint that to our knowledge has not yet been 
discussed. Let us define the Fourier decomposition as 

m(t) = m + acos(ojt + <p{) + bcos(2ujt + fo) + • • • 

and as usual 4>2i = 02 — 20i. 

The lightcurves over a whole Blazhko cycle have been published by Walraven 
(1949) for RR Lyr and by Balazs & Detre (1939) for AR Her. In Fig. 5 we show 
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Figure 6. Schematically, location of strange Cepheids in HR (left) 
and P-L (right) diagrams. 

the variation of the pairs of Fourier parameters, a, b, and foi- All quantities are 
seen to oscillate about some center, which is of particular interest for the phase 
4>2i- The fact that the phase does not run through 2ir over a cycle imposes a 
severe constraint that can be used to eliminate some models. 

5.3. The dip in the Galactic Cepheid period histogram 

In 1977 Becker, Iben & Tuggle published Cepheid period histograms for several 
galaxies. The histogram for the Galaxy and for M31 showed a pronounced dip in 
the 8-10 day period range, whereas the corresponding histograms for the LMC 
and SMC were devoid of such a deficiency. Trying to explain the dip on the 
basis of their stellar evolution calculations Becker et al. had to invoke an ad 
hoc double-humped birthrate function. Nonlinear calculations show that this 
is no longer necessary (Buchler, Goupil & Piciullo 1997). Indeed, a perhaps 
unexpected side effect of the new opacities is that the fundamental limit cycle 
of the Cepheid variables can be unstable. This instability is found to occur in 
the 8-10 day period range for metallicity parameters 0.013 < Z < 0.035. Note 
that this is consistent with a dip in the Galaxy and M31 and with the absence 
of a dip in LMC and SMC. 

5.4. Strange Cepheids and RR Lyrae 

It has recently been found that strange modes can occur even in weakly nona- 
diabatic stars such as Cepheids and RR Lyrae. A thorough study of the phe- 
nomenon has shown that these modes are surface modes that can be self-excited 
to the hot side of the blue edge of the normal Cepheid instability strip (Buch- 
ler, Yecko & Kollath 1997). The strange modes are recurrent at higher order 
wave- vectors, but the lowest ones have typical periods 1/4 to 1/5 that of the 
fundamental pulsational mode, i.e. they have periods ranging from ~0.2 days to 
~10 days, depending in their luminosity. Their locations in schematic HR and 
PL diagrams are shown in Fig. 6. 

What does one expect the pulsations to look like? We have computed some 
nonlinear (radiative) models and find limit cycles with amplitudes in the mil- 
limag and 10-100 m/s ranges, respectively. It might be feared that, because the 
strange modes are surface modes their driving could be destroyed by convection. 
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Preliminary computations with the turbulent diffusion hydro code however in- 
dicate otherwise. 

Finally we note that the same trapping and driving mechanisms also work 
in RR Lyrae models and that therefore strange RR Lyrae should also exist on 
the hot side of the RR Lyrae instability range. 

6. Conclusions 

The theoretical study of stellar pulsations is still faced with many challenges. We 
have seen that radiative hydrocodes, while giving decent agreement with many 
observations, are not fully satisfactory. We hope that a proper inclusion of the 
important dissipative effects of turbulent convection will help resolve many of 
the extant difficulties and discrepancies. 
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